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I. INTRODUCTION 



In previous work, effective field theory (EFT) [0, |, |, |, H, @ was applied to the normal 
ground state of a dilute Fermi system with short-range interactions 0,11] and subsequently 
extended to finite systems as a realization of Kohn-Sham density functional theory (DFT) 
@, EJ These studies are part of a program to go beyond mean-field approaches to finite 
many-body systems while maintaining their computational simplicity and phenomenological 
strengths. Our ultimate goal is to systematically describe medium and heavy nuclei using 
EFT and renormalization group methods, including nuclei far from stability where EFT 
power counting and error estimates are desirable to enable reliable extrapolations. For nuclei 
as well as trapped systems of atomic gases, attractive interactions cause the normal ground 
state to be unstable to the formation of Cooper pairs. In the present work, we consider an 
extension of the EFT formalism for uniform systems to allow for pairing correlations in a 
way that is also consistent with Kohn-Sham DFT. 



Our strategy is to apply an effective action formalism [12], [13|, [14| to calculate the free 
energy order by order in an EFT expansion. The effective action formalism is often used 
in condensed matter applications to discuss superconductivity, star ting from a path inte- 



gral with a four-fermion contact interaction (see, for example, [la . |l6|). In the standard 
approach, a pairing field is introduced as a charged scalar auxiliary field and shifted to elim- 
inate the four-fermion term, leaving only a Gaussian integral over the fermion fields. After 
performing this integral, a conventional effective action is derived by a functional Legendre 
transformation with respect to the auxiliary field. A minimum at a nonzero expectation 
value of this field (which is proportional to the pairing gap for a short-range interaction at 
leading order) indicates spontaneous symmetry breaking of the phase symmetry related to 
fermion number conservation, which implies the normal ground state is unstable to pairing 
and determines the "superconducting" ground state. 

Instead of introducing an auxiliary pairing field, we couple an external source (or, more 
generally, a pair of complex conjugate sources) to the local, composite pair density [Eq. ([19]) ]. 
In general this source would be space-time dependent, but for a uniform system in its ground 
state, which we consider here, the source can be treated as a constant. This external source 
is analogous to the source introduced for finite systems in Ref. 0, which was coupled 
to the composite fermion density. As in that case, we carry out a functional Legendre 
transformation order by order in the EFT expansion via the inversion method introduced 
by Fukuda et al. 13, 11] (see Ref. [l9| for an application to superconductivity). If we keep the 
space dependence of the sources and the corresponding densities, we obtain a generalization 
of the EFT-based density functional theory (DFT) formalism of Ref. [qJ (see also Refs. [2(J 



2ll . |22| . l23l . [2J, |25j). In the present work, the limitation to uniform systems means that the 



inversion method parallels the construction used long ago by Kohn, Luttinger, and Ward 
2^, 27] to show the relationship of zero-temperature diagrammatic calculations to ones using 



the finite-temperature Matsubara formalism in the zero-temperature limit. 

In Refs. 28] and [2^], density functional theory for superconductivity was proposed using 
a nonlocal source coupled to the nonlocal pair density. While this approach can be embedded 
in an effective action framework, we first want to explore the use of local sources, which are 
greatly preferred for finite nuclei from a computational point of view. In particular, we seek 
a consistent generalization of the Skyrme-Hartree-Fock-Bogoliubov approach (30| . 

Potential problems arise, however, because the use of local composite operators and zero- 
range interactions leads to new ultraviolet divergences even at the mean-field (Hartree-Fock) 



2 



level when pairing is included. Divergences at this order are not encountered when coupling 
an external source to the fermion density V V> ^ut a PP ear now because the composite 



operators ipip and ip^ip^ need additional renormalization 3JJ. The divergences at leading 
order are symptomatic of generic problems identified long ago by Banks and Raby [32| that 
arise with effective potentials of local composite operators. (Such divergences would be 



avoided by working with nonlocal sources as in Refs. [28l . |29I|.) These problems inhibited 
for many years the use of effective actions of composite operators although more recently 

M 



Verschelde et al. 33, 34, 35 



36|, [3j 



have revived their use for 



and Miransky et al. 

relativistic field theories. In this work, we show how to remove the new UV divergences 
through next-to-leading-order (NLO) by adding a term quadratic in the pairing source, 
using a renormalization prescription compatible with the inversion method and conventional 
treatments of pairing at leading order. 

The BCS gap equation for contact interactions has conventionally been renormalized 
at leading order (LO) by using the equation for free-space scattering, which has the same 
divergence [see Eqs. f l5^|) and f[6"Uj) ]. The four-fermion coupling is eliminated in favor of 
the scattering length, leaving finite, renormalized equations. This approach was applied in 
Refs. 39[] an d |40l| to obtain analytic formulas for the pairing gap and ground state energy at 
LO. In Ref. [40|], dimensional regularization with minimal subtraction (DR/MS) was used. 
DR/MS was also used to good advantage for the dilute Fermi calculations of the normal 
ground state in Ref. 0]. In extending the formalism of Ref. JtJ to include pairing, we use 
a more general subtraction scheme (called "power divergence subtraction," or PDS), that 
lets us verify explicitly that the generating functional, effective action, gap equation, and 
observables are renormalized order by order in the inversion approach. We give renormalized 
expressions up to next-to- leading order (NLO), which contain the universal correction to the 
pre-factor of the gap at low density, first calculated in Ref. (see also Ref. j4~2|). 

The plan of the paper is as follows. In Sect. [Til we briefly rederive the EFT for a dilute 
Fermi gas in the effective action formalism using the inversion method of Kohn, Luttinger, 
and Ward [26}. In Sect. IHH the approach is generalized to include pairing by introducing 
a source coupled to the pairing density and then carrying out the Legendre transformation 
to get an effective action as a functional of fermion and pair densities. The construction of 
the EFT effective action functional via the inversion method starts with the zeroth order 
functional. In Sect. HVl we carry out the leading order (LO) calculation in detail, reproducing 
standard results and extending the dimensional regularization results of Ref. [40(]. In Sect. |V] 
we extend the calculation to next-to-leading order (NLO), which introduces the so-called 
"induced interaction" . Section |VT] is a summary and some additional details are included in 
the Appendices. 



II. INVERSION METHOD WITHOUT PAIRING 
A. Lagrangian and Generating Functional 

In contrast to the zero-temperature, Minkowski-space EFT treatment of dilute systems in 
Ref. [3] , we work in Euclidean space in the Matsubara (imaginary time) formalism at inverse 
temperature f3. We use a four-vector notation x = {r, x}, with the r dependence usually 
implicit, or k = {a; n ,k} in momentum space, where u n is a fermion Matsubara frequency. 
Throughout this work, however, we will consider the zero temperature (/3 — > oo) limit, for 
which the differences from Ref. [3] are minor (and the results equivalent after frequency 
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integrals are performed). 

The EFT Lagrangian for a nonrelativistic spin- 1/2 fermion field with spin- independent 
interactions is the same as in Ref. [7J], but converted to Euclidean form (which means t — > —it 
and an overall minus sign in defining Ce and the interaction Lagrangian C 



d_ 



2M. 



$ + ^° (^V) 2 - — IW^Vfy) + h.c. 



16 



d_ 



2M 



(1) 



where V = is the Galilean invariant derivative and h.c. denotes the Hermitian 

conjugate. 1 The terms proportional to C 2 and C 2 contribute to s-wave and p-wave scattering, 
respectively, while the dots represent terms with more derivatives and/or more fields, as well 
as renormalization counterterms. The Lagrangian Eq. ([T]) is a particular canonical form, 
which can be reached via field redefinitions. For example, higher-order terms with time 
derivatives are omitted, as they can be eliminated in favor of terms with spatial derivatives 
(see, for example, Ref. j43[). We will restrict ourselves here to spin-1/2 (i.e., degeneracy 
v = 2) and spin-independent interactions, which is sufficient to illustrate the formalism 
and the renormalization issues; generalizing to isospin and spin- dependent interations is 
straightforward. We also consider only the Co vertex here, but comment at the end on the 
inclusion of vertices with derivatives. 

We first consider a Euclidean generating functional with chemical potential /i and external 



sources r\{x) and iy{x) [16l. l44j: 



Z[ri,rf;n] = e 



-W[rirf ;n] 



t p - f d4x [Ce - fiipi(x)rp a (x)+rji(x)^ a (x)+^l l (x)r] a (x)] 



(2) 



where J d A x includes a dr integration that runs from — {3/2 to {3/2 (to facilitate the (3 — > oo 
limit) and anti-periodic boundary conditions are imposed. For simplicity, normalization 



factors are considered to be implicit in the functional integration measure. (See Refs. [17, 



18l . I44j for more detailed treatments of similar path integrals.) We have written the spin 
index a explicitly in Eq. (|2J); we will interchangeably use a = {1,2} and a = {f, 1} in the 
sequel to denote individual spins, and the spin indices will be left implicit where there is no 
chance of confusion. We have kept the chemical potential \x separate from the Lagrangian 
in Eq. ([2]) to emphasize its role as an external source. In the next section, we will add a 
corresponding source coupled to the pair density. 

A conventional perturbative expansion is realized by removing the interaction terms from 
the path integral in ([T]) in favor of functional derivatives with respect to r\ and and 
performing the remaining Gaussian integration over ip and ip^ (l6l . 44j: 



Z[r},rf;fj\ = Z e 



(3) 



where the spin indices are implicit and we have introduced the noninteracting partition 
function Z . Explicit expressions for the Matusbara Green's function Q in coordinate and 



1 We use units with h = 1 . 
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FIG. 1: Hugenholtz diagrams for the unrenormalized Wi functionals in a homogeneous, dilute 
Fermi system. The vertices are defined in Fig. [5j 



momentum space can be found in Ref. [44|. The linked-cluster theorem [44j shows that the 
difference of the interacting and noninteracting thermodynamic potentials Q and Qq is given 
by the sum of connected diagrams from the expansion of Z, with the external sources 
and rj set to zero at the end: 



T, V) - Q Q (fj„ T, V) = ±(W[0, 0; fi] - W [0, 0; fi)) 



(4) 



The connected Feynman diagrams that sum (with appropriate symmetry factors) to W — Wq, 
which are labled Wi with % > 1, are shown in Fig. [H w ith the fermion lines representing 
Matsubara propagators Go with chemical potential /i [44] . The subscript i labels the leading 
order to which Wi contributes in the dilute EFT expansion. (In the DR/MS scheme, each 
diagram contributes to only one order but this is no longer true in DR/PDS.) The chemical 
potential can be determined for given N by inverting the thermodynamic expression for the 
particle number, N = dQ/dfi. The regularization and renormalization of divergences arising 
in the evaluation of the Wi are described below. 



B. Kohn-Luttinger-Ward Inversion 



The Kohn-Luttinger-Ward (KLW) theorem 26|, 27] relates the perturbative calculation of 
diagrams using the finite-temperature Matsubara formalism in the zero-temperature limit to 
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the calculation of diagrams using zero-temperature perturbation theory, which was applied to 
the dilute Fermi gas in Ref. Q. In particular, the diagrammatic expansion of f2(/i, V, T — > 
0) is, under specified conditions, the same as the free-energy F(N, V,T = 0) evaluated 
using zero-temperature propagators with the non-interacting chemical potential /i and with 
anomalous diagrams omitted (see Ref. [?J for Feynman rules). The anomalous diagrams 
are those such as Fig. [T](c) and (d), which vanish in zero-temperature perturbation theory 
because they have particles and holes in the same intermediate single-particle state (i.e., 
factors of 9(fi — e a )9(e a — /iq) = appear). The conditions of the theorem require that the 



symmetry of the non-interacting and interacting systems agree [4J, |45[ . 

A demonstration of the KLW theorem using an inversion method for the case of an 
electron gas is presented in Ref. (45[. We adapt this demonstration to the case of an EFT 
with a well-defined expansion. While we illustrate the procedure with the expansion relevant 
to a natural short-distance interaction (which is an expansion in the Fermi momentum k-p 
times the effective range parameters a s , r s , and so on [3]), the discussion is more general, 
requiring only a hierarchical expansion. Thus, for example, a nonperturbative (in diagrams) 
1/N expansion works just as well. 

The basic plan is to carry out the conventional thermodynamic Legendre transformation: 

F(N) = n(fj) + fj,N , (5) 

with /i(iV) obtained by inverting N(fi) = — {dVL/d^TV-, which determines the mean- number 
N of particles. We expand each of the quantities in Eq. ([5]) about the non-interacting system: 

Ufa) = n (//) + + + • • • , (6) 

fj. = fa + Hi + fi 2 H , (7) 

F(N) = F Q (N) + Fi(JV) + F 2 (N) + • • • , (8) 

where the subscript indicates the order of the expansion. (Note that the subscript is just 
a counting parameter that does not have to correspond to a power series in the expansion 
parameter; e.g., in the Coulomb case the expansion parameter is e 2 but Q 2 has both an 
e 4 term and the correlation energy of order e 4 lne.) The non-interacting system refers to 
a system of zeroth order in the EFT expansion parameter. This means the zeroth-order 
system has no internal interactions among the particles, but it can include external sources 
(we exploit this freedom below). Each is given by (l//3)Wj[0, 0; //], where the Wj's 

correspond to the diagrams in Fig. [TJ 
The number of particles is 

N(uTV) = -(—} =-^2-^--^ + ... (9) 

Note that we could simply use the unexpanded first equality in Eq. (jUj) together with the 
series in Eq. since they define a parametric relation between N and Q in terms of /i 



451 ] . Since we have in mind the extension to finite systems, we pursue an alternative that 
generalizes to Kohn-Sham density functional theory (DFT). To carry this out, we perform 
the inversion in Eq. (jH]) order by order in the EFT expansion, treating iV as order zero in 
the expansion. (That is, we ensure there are no corrections to iV at higher order.) This 
means that the zeroth order equation, 

[ttI • ^ 
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is the only equation to which N contributes (by construction). Thus the "exact" N is ob- 
tained from the non-interacting thermodynamic potential. This might not sound impressive, 
but the analogous situation holds when we generalize \i to be position dependent or coupled 
to the pair density in a finite system. In these cases, it is the exact, spatially dependent 
fermion or pair density (with appropriate renormalization conditions, see below) that is 
obtained from a non-interacting system with a single-particle potential that is the general- 
ization of fj,Q. This is precisely the description of the Kohn-Sham system (see Refs. |9j and 
(To| for details on carrying out the DFT case without pairing). 

Equation ([TO]) determines N(fio) at any temperature, from which we can find fJ>o(N) for 
any system for which we can identify Qq (the inversion is unique since fio is a monotonic 
function of N j45[). If we have a uniform system with no external sources, /io is the chemical 
potential of a noninteracting Fermi gas at temperature T with density N/V. In particular, 
at T = with no external potential and spin-isospin degeneracy u, 



^(N) = (6n 2 N/vV) 2/3 = kl/2M = e° F . 



(11) 



The first-order equation extracted from Eq. (j5J) has two terms, which lets us solve for \i\ 
in terms of known (from diagrams) functions of 







djj, 



+ Hi 



d\i 2 



/'i 



[dQx/dfj] 



[d 2 n /dfi 2 ' 



(12) 



At second order, we can isolate and solve for fi 2 , eliminating \i\ using 012]) . This pattern 
continues to all orders: /ij is determined by functions of /x only. 
Now we apply the inversion to F = Q + jiN: 



F{N) = n (n )+fx N + /n 



F 



+ QtQio) + HiN 



=/*() 



Fi 



+ A«2 



1 2 



d 2 Q c 
dfi 2 



+ Q 2 M+fM 2 N + .-.(13) 



-v— 

Fi 



The HiN term always cancels with ^i[dflo/dfi\^ =lJi0 in Fi forz > 1 because of Eq. (fTOl) . leaving 

+ ••-, (14) 



1 [9fii/^ =M0 



F(N) = F (N) +S M + n 2 (, o) 2[d2W 

Fi v 



/•»=A»0 



F 2 



where we've also used Eq. 0121) to simplify F 2 . The expansion for F can be extended 
systematically, but this is all we need here. (Higher orders can be found by following the 
prescription in Refs. [H, 47].) 

This construction is rather general. The Kohn-Luttinger-Ward theorem explores a par- 
ticular case, the T — > limit, in which the second term in F 2 in Eq. ffT4"l) cancels precisely 
against the anomalous diagram in f^, as illustrated in Fig. [2j This cancellation of deriva- 
tive terms and anomalous diagrams occurs to all orders in the expansion. (An analogous 
cancellation was noted in the Kohn-Sham density functional (DFT) inversion in Ref. [9| 
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oOo+oCK>o 




FIG. 2: Cancellation of the anomalous diagram at NLO. The double lines represents the inverse 
of [<9 2 fV<V]» in Eq. (HU). 

and appears again in Sec. IIII A[ ) The end result is an expression for the free-energy F(N) 
in terms of the diagrams used for fli(p), only evaluated with p = p and excluding the 
anomalous diagrams (both of which simplify the evaluation of F(N)\). This is precisely the 
formalism used in Ref. [7| for a uniform low-density Fermi gas at zero temperature, where 
p appeared as the Fermi energy of Eq. ffTTl) . 

If we think of p as an external (constant) source coupled to the fermion number, we 
can also imagine having sources coupled to other quantities and making analogous inver- 
sions. The generalization to multiple chemical potentials (i.e., constant sources coupled to 
conserved charges) is straightforward but adds nothing new. There are, however, two more 
interesting generalizations: 

1. Add a second "chemical potential" j coupled to a non-conserved charge, with j equal 
to zero in the ground state. The same inversion method holds, but since we don't 
have a conserved charge such as N, we can't determine jo from the analog to Eq. ffTUl) . 
Instead, we use jo + ji + ji + ■ ■ • = in the ground state to solve for j iteratively. (That 
is, guess j , calculate ji for i > in terms of j , find a new j from j = — (ji + J2 + *••), 
and repeat to self-consistency.) It is this method that we use here to incorporate pairing 
for a uniform system with no external fields. 

2. Add a spatially dependent (classical) source J(x) coupled to the density operator ^i/}, 
which has expectation value p(x) = (ip^ip). Then construct the functional Legendre 
transformation F[p] = Q[J] — f pj, where J(x) = <5F[p]/5p(x) — >• in the ground state. 
This leads to the Kohn-Sham DFT discussed in Ref. [9|, where F[p] is proportional to 
the effective action T[p] for the composite operator ip^ip. 

In Ref. [Tcj, a spatially dependent source coupled to the kinetic energy density is included; 
in this case the inversion method leads to a generalization of conventional Kohn-Sham DFT 
that is analogous to the Skyrme-Hartree-Fock energy functional [48]. A natural generaliza- 
tion of Skyrme-Hartree-Fock-Bogoliubov will follow by extending the present discussion to 
J(x) in finite systems. 



C. Renormalization 

Starting with W^l/t], the diagrams in Fig. [T] may have ultraviolet divergences. In Ref. Q, 
these divergences were regulated and renormalized using dimensional regularization with 
minimal subtraction (DR/MS). This choice provided a clean factorization of each diagram, 
with the dependence on the effective range parameters solely in the coefficients and the 
integrals depending on the Fermi momentum kp times universal geometric factors. The 
result is a very systematic power counting, with each diagram contributing to exactly one 
order in the EFT expansion. 
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However, while DR/MS is convenient and efficient and correct, it can obscure the renor- 
malization process. In the present case, since there are subtle questions about renormalizing 
effective actions with composite operators, we adopt the power divergence subtraction (PDS) 
scheme of Kaplan, Savage, and Wise which is a generalization of DR/MS that was intro- 
duced to provide a consistent power counting in the two-body problem when the scattering 
length a s is unnaturally large 0, 5(|. (An alternative for that problem is a momentum 
subtraction scheme proposed by Gegelia |5lfl.) A parameter is introduced, which we call A 
to avoid confusion with the chemical potential (this same parameter is called \i in Ref. 49|). 
as part of the finite subtraction. To apply PDS for natural scattering lengths, we take 
A <C l/a s and expand consistently for small A. The DR/MS scheme is recovered by taking 
A = 0. The point of keeping A explicit is that we can verify that we have successfully 
renormalized (to a given order) by observing the cancellation of A dependence between the 
couplings and loop integrals. 

For example, consider Fig. [Tfb) evaluated with fj, . After applying zero temperature 
Feynman rules, carrying out the frequency integrals, and transforming to center-of-mass 
variables, we obtain (see Ref. 0] for further details): 



£ 2 = AClMv{v - 1 



(2tt) 3 
x [1 - 6(k F - 



dH 



3-D 



U 



\)][i-e(h 



d D u 
(27p 



6(h 



t\)6{h 



s — u 



t 2 — u 2 + ie 



|s-t|) 

(15) 



with v — 2 the spin degeneracy. The integral over u in Eq. (j!5j) . which is regulated in D 
spatial dimensions, contains a linearly divergent term for D = 3 coming from the piece of 
the integrand without any 6(k F — |s±u|) factors. In PDS, with the conventions of Ref. 49 
this integral is: 



3-D 



d 



D 



U 



PDS 



u 2 + ie 



An 



(A + it) 



(16) 



(2n) D t 2 

The imaginary part cancels with an equal and opposite imaginary part from the remainder 
of the (finite) integral over u [identified using l/(t 2 — u 2 + ie) = V{l/(t 2 — u 2 )} — in5(t 2 — u 2 )], 
so that the energy density is real. The real parts of the integrals that contain 6(kp — |s± u|) 
factors reproduce the result for £ 2 given in Ref. Q. The contribution proportional to A is 
easily evaluated after changing variables to q = s + t and q' = s — t: 



r d 3 <* r (fit 1 
5S 2 = -AClMu{u-l) J-^-e(k F -\s + t\)9(k F -\ S -t\)-A 



47T 



--L C 2 M(l-l/v)A 

S7T 



V 



d 3 q 



-, 2 



--C 2 M(1 - l/v)Ap* 

07T 



(17) 



where p is the fermion density. Thus, we have an apparent dependence in the energy on the 
arbitrary parameter A. 

However, we also have a contribution at 0(A) from the diagram in Fig. QJa) from the 
dependence of Cq on A in the PDS scheme. In particular, the requirement of A independence 
for the free-space scattering problem yields |49| : 



C (A) 



47T 

M 



1 



-A + 1 /a s 



Ana, 



Ana 2 . . n . . 2 
+ -jf^ + 0(A 2 



cp + cr + ■ ■ 



(2) 



(18) 
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(2) 

FIG. 3: Renormalization of the bubble in free space determines Cq to remove the linear A 
dependence. This implies the renormalization of the beachball (or ABA basketball 52() diagram 
Fig.QTb) at 0(A). 



W[fi,j] 



+ 



x 

3 




----x 

3 



+ 



FIG. 4: Feynman diagram at second order in a perturbative expansion of W\p, j] in j(x) 



If we substitute the leading (independent of A) term in Eq. ffl~8l) for C in Eq. f|T7|) . we get 
a term proportional to kpa 2 s A. This cancels precisely with the contribution for the Hartree- 
Fock diagram Fig. CQa) evaluated with a Cq vertex (as implied by Fig. [3]), which removes 
the linear A dependence. Higher-order insertions of Co cancel between these diagrams and 
higher-order diagrams in the k-p expansion. 



III. EFFECTIVE ACTION AND INVERSION METHOD WITH PAIRING 

A. Generalized Generating Functional and Inversion 

To allow the probing of pairing correlations, we extend our generating functional to 
include the source j(x) coupled to the pair density, 



In general we should include both j(x) and j*(x), but in the present case it suffices to 
consider real and constant j. The Grassmann sources rj and are not shown, but they can 
be added to generate a perturbative (i.e., diagrammatic) expansion for Z or W. 

As discussed in Ref. [53| (see also Verschelde et al. 33, 34, 3^|), to make finite a gener- 



ating functional with a source coupled to a local composite operator, we will need to add 
as counterterms the most general linear combination of vertices of non-positive dimension 
made from that operator and the source. A perturbative expansion of W in powers of j 
immediately reveals a divergence, shown diagrammatically in Fig. HI To renormalize in the 
present case, we include the term |0 2 ; which plays a role in the renormalization (see below) 



analogous to the quadratic charged scalar term in the auxiliary field approach [15j,ll6|. We 
note that this term is also the source of potential problems with the interpretation of the 
effective action [32 ]; we comment further on this issue in Sect. IIVI Other counterterms in 
Eq. f|T9|) . including one proportional to + V'J.V't)' are implicit. 
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The generating functional implies that for constant j and a uniform system, 

4, = + + 0') M - = ^^j 1 ^ • (2i) 

Note that corresponds to the usual definition in the ground state, where j = 0. The 
explicit dependence of on j through ( is a consequence of the need to renormalize 



oep 

33l . l34j . |35I |. We generalize the KLW derivation from Sect. Ill Bl by carrying out a double 



Legendre transformation, making the associations 



F^^T and Q^^W, (22) 

P P 



to switch variables to p and 



^r[p,0] = -±^w\n,j]+np-j4> . (23) 

We proceed in analogy to Eqs. (J6j) (JHJ) to construct the effective action T[p, 0] using the 
inversion method. In particular, we now have the expansions 

W\fiJ] = W [p,j] + W 1 [p,j] + W 2 [p,j] + • • • , (24) 

p = p + pi + p 2 H , (25) 

J = jo+h + 32 + ■■■ , (26) 

T[p,j] = T [p,j]+T 1 [p,j] + T 2 [p,j} + --- . (27) 



As before, the inversion is carried out by matching Eqs. (1201) . ( 1211) and ( 1231) order by 
order using (|24l) - (l271) . with both p and <fi counted formally as zeroth order. A complication 
in carrying this out is that will in general be renormalized order by order, which implies 
contributions to the left side of Eq. (T2TT) beyond zeroth order. To avoid this, we will use the 
freedom in the renormalization conditions to require that there are no corrections to 0, i.e., 
that and p are given exactly by the zeroth-order equations, 

1 dW [po,jo] 1 dW [po,j ] 

9 ~ ~W dp ' ^'W d 30 ' (28) 

where it is implied that po and jo are determined as functions of p and 0. The contributions 
to Ti from p^p and jj0 when % > 1 cancel as in Eqs. f[T3"j) and CHI) after using Eq. (T2"g|) . The 
anomalous diagrams in W 2 cancel as in Eq. ( 1T41) and Fig. [2] against additional terms from 
the inversion (§1, 1(|, leaving (with W 2 denoting the non- anomalous part) 

T[p,<j>) ^T [p, ( f>}+W 1 [p ,jo} + W2[po,3o}+--- ■ (29) 

r\ r 2 

Thus, in practice we only calculate the non-anomalous diagrams to this order. 

To carry out a calculation, we truncate the expansion at a given order, construct the Wi 
needed in the corresponding truncated version of Eq. ( 1291) . and determine the pi and ji for 
i > 1 from 

1 <9r,[ P ,0] . _ _ 1 gr t [ P ,0] 

dp ' 3l ~ PV d<P ■ {M) 



11 



+ V + V + 



P/2 + k P/2 + k' 

X - 

P/2-k P/2-k' 

-<k'|V EFT |k) -Co -C a - ' - -t^k-k' 






k 2 + k' 2 




= y\ + - 

--Do 

FIG. 5: Feynman rules for — (k' | V^jft I k) and the leading 3-body contact interaction. The spin 
indices have been suppressed. 

The values po and jo are fixed by requiring both that a specified p to be reproduced by 
Eq. (128p and that j = 0, which means 



Jo = -(ji+j2 + ---) • (31) 
In general, these conditions can only be established numerically. 

B. W [p,j] and T \p,<f>] 

To find Wo and To, we need to solve a non-interacting system in the presence of external 
sources po and jo- The path integral expression for VFof/^o, jo] with constant sources is 

Z [p J ] = e - W °^ujo] = f DlpaD^l e -/ rf4 ^[^(^r-V 2 /2Af-Mo)^+joW T ^+^^|)+K (0) Jo] ; 



(32) 

where we restrict the a sum to run over spin up and spin down and the x dependence of 
the ip fields is implicit. We denote the value of £ at this order as which should depend 
on A only. Since the path integral is quadratic in the fermion fields, it can be evaluated 



exactly in the Nambu-Gorkov formalism as Trh^Gj^ 1 ) (see, for example, Refs. [la . Il6|). 
where G^, 1 is the quadratic part of the Lagrangian and Gk s is the associated Kohn-Sham 
Green's function. Alternatively, a canonical Bogoliubov transformation solves the problem 
exactly and is particularly transparent and familiar. (In a finite system, we would solve 
the zeroth-order system by solving Hartree-Fock-Bogoliubov equations with external Kohn- 
Sham sources corresponding to a spatially dependent chemical potential and a spatially 
dependent pair density.) This is carried out in Appendix [A] and we simply quote that result 
for Wo: 

1 f rl 3 k 1 

-pjWoW Jo] = Jj^ys^- E k ) + ~ 2 C m jl , (33) 

where 

^ = ei-po, (34) 
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with = k 2 /2M and 



E k = ^e k +3l- (35) 
The zeroth-order effective action is then 

1 1 

Equations for p and follow from Wo: 

N 1 dW f d 3 k ( i k 



and 



P ~V~ f3Vdp 'J (2vr)3 y 1 E J ' (37) 



1 dW f d 3 k jo = 



+ C w io = 0b + C W io • (38) 

PV dj J (2tt) 3 E fc 

Except for the extra term proportional to these results are the same as from standard 
BCS calculations, but with jo appearing in the place of a constant gap. This is illustrated 
further in the left panel of Fig. [HI which plots quantities defined in Eq. flA14[) . With those 
standard definitions, p = J d 3 kvl and 0b = <ps / d 3 kukVk- 

The integral on the right side of Eq. (138]) . which we have defined as the "bare" pairing 
density 0b, is linearly divergent (see the right panel of Fig. [9]). This requires regularization 
and then renormalization through an appropriate choice of To apply the DR/PDS 

scheme, we have adopted and extended the integrals derived for DR/MS in Ref. 4(j in 
Appendix [Bl From Eq. (1B6I) we have immediately that 

= -jo/(O)+C (O) Jo • (39) 

Since the interaction (i.e., C ) should not appear explicitly in this zeroth-order equation, 
requiring to be the full renormalized result independent of A uniquely determines 

C<°> = ^ (40) 

by dimensional analysis. This results in an analytic expression for in terms of Legendre 
functions [see Eqs. (1B5I) and (1B6I) ]: 

= -jo 7(0), (41) 

which is manifestly independent of A. Upon applying ( ]B6j) term-by-term to the density 
equation, we find (as expected) that the A dependence cancels: 

P = (|) J (l - = - 1(1) + a*o 1(0) = -HI) + Mo 7(0) . (42) 

By applying Eqs. (IB 71) . these expressions for and p can be expanded at small x = jo/po 
[see Sect. ED and Eq. flBTl) ]. 
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C. Feynman Rules for Wj 



The Feynman rules for the contributions to Wi[/j,,j]//3V = 0,%/V from the diagrams in 
Fig. [D are similar to those given in Ref. 0] (which in turn take the standard form, e.g., see 
Ref. 45]), but with Nambu-Gorkov matrix propagators. These propagators arise from the 
diagonalization of Wq. The only difference from standard treatments of superconductivity 
(e.g., see is that j appears instead of a constant gap. 

The basic Feynman rules in the zero-temperature limit are as follows. First, draw all 
distinct, fully connected diagrams contributing to a given order in the EFT expansion of 
W. (In this paper we determine diagrams according to the k-p expansion for a dilute gas 
without pairing, so the power counting is the same as in Ref. [3].) Distinct diagrams are 
those that cannot be deformed to coincide with each other, including the direction of arrows. 
To evaluate a diagram: 

1. Assign nonrelativistic four-momenta (frequency k and three- momentum k) to all lines 
and enforce four-momentum conservation at each vertex. At finite temperature, the 
frequency is a discrete Masubara fermion frequency [IF 



For each vertex, include the corresponding term from the two-body effective potential 
— (k' J Veft I k) (or analogous three-body interaction), as shown in Fig. [5J (Recall that 
k and k' are relative momenta.) The Nambu- index structure of the vertices has been 
suppressed in Fig. [51 For spin-independent interactions, the two-body vertices have 
the structure [(r 3 ) a -y(T 3 )i3s ± (r^aslrs)^}, where a, (3 are the Nambu indices of the 
incoming lines and 7, 5 are the Nambu indices of the outgoing lines. The plus sign 
applies to Cq and C2 vertices and the minus sign to C 2 vertices. 

For each internal line include the matrix propagator — Go(fc), where k = (ko, k) is the 
four-momentum assigned to the line, and 



Go(fc) 



£o(fco,k) J" (^o,k) 
J$(k ,k) -Q (-k Q ,k) 



with 



Go(k) 

and 

r Q (k) = rt{k)-- 



4 



+ 



1 



ik - E k ik + E k 2E k 



E k + £fc E k — £ k 



-u k v k 



1 



1 



ik - E k ik + E k 



ik - E k ik + E k 



Jo 



(43) 



(44) 



2E k 



ik - E k ik c 



Ei. 



(45) 

(Note: These conventions agree with those in Refs. [45j and [54| and differ by a sign 
from those in Ref. Hj].) The u k s and v k s are defined in the Appendix in Eq. flA14j) 
and are seen here to correspond to the residues at the poles of the noninteracting 
Green's functions. 

3. Perform summations on Nambu indices in the diagram (e.g., do the matrix multipli- 
cation). Include —1 for every closed fermion loop (i.e., for each trace). 

4. Integrate over all independent momenta with a factor J d A k/(27i) 4 where d A k = dk d 3 k. 
[At finite temperature, the integral J dk /2n is replaced by the sum l//3J2 n over 
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fermion frequencies u n [451].] If a spatial integral is divergent, it is defined in D 
spatial dimensions and renormalized using the PDS subtraction scheme as discussed in 
Ref. [49] . For lines originating and ending at the same vertex, multiply by exp(±ifco^) 
and take the limit rj — > + after the contour integrals have been carried out. Use 
exp(ik rj) for Q ^ and Tqj. and exp(— ik rj) for Q Q ^. This procedure ensures that the 
correct term in the Green's function is picked up. 

5. Multiply by a symmetry factor — l/($n!=2 X (^) m ) wnere S is the number of vertex 
permutations that transform the diagram into itself, and m is the number of equivalent 
/-tuples of lines. Equivalent lines are lines that begin and end at the same vertices 
with the same direction of arrows. 



IV. RENORMALIZATION TO LEADING ORDER 

In this section, we carry out the inversion procedure to leading order (LO), which means 
calculating Ti[p, <f)] and determining equations for fio and jo- We will work here in the zero- 
temperature limit, with frequency integrals rather than frequency sums. The final results 
are well known, although to our knowledge this precise formalism has not been used. 

Since i\[p, <p] = Wi[fi [p], j [4>}] (see Ref. [§]), we first find Wi[fi ,jo] by applying the Feyn- 
man rules to the diagram in Fig. QJa). The overall symmetry factor is 1/2 (one equivalent 
pair of lines), so the expression is (after cancelling signs), 

gyWxifio, jo] = [{T 3 ) ay (T 3 ) ps + {r 3 ) aS (T 3 )^} 

X /(0/|^ G ^ G ^ 

1 (1) f d*k [d*k> ~ ~ 

~ 2 J (2tt)4 J (2~7r)4^ * k ' ^ k ' 

= Cf ) / 1^ [Qa&r v + *3s*bp] + ^C (1) J 2 , (46) 

where we've suppressed the e ±tkoV convergence factors and (except for i^ 1 )) additional coun- 
terterm corrections discussed below. The individual integrals we need are: 

d \ 5o ( ± *o, k) ^ = = [***Llk = l p , (47, 



(2vr) 4 uv u ' ' J (2tt) 3 k J (2tt) 3 2E k 2 

|L^,, k)c ,.„ . ||L^._ J0 (A) 3 - D y |i_i_4, B , (48 ) 

(and JF = JF^). The integral in Eq. (14"51) has a linear divergence, which we've regulated in 
D spatial dimensions (and identified 0b accordingly). 

In addition to the |C^Jo counterterm in Eq. (146 p . there is a counterterm SZ^jtyfa^ + 
iplipl), which contributes to W\ a term proportional to Jq4>b (see Fig. [6]). Our renormalization 
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OO + x--"0 + 0"--x + X— D--X 



a 



(i) 



jo $Z 



(i) 



M 1} Jo 



Jo 



C (1) Jo 



FIG. 6: Diagrammatic representation of the renormalization of T\ 



prescription is to choose this counterterm and so that we convert 0b to the renormalized 
from Eq. fl38|) [with £W from Eq. (j4"U|) ]. This takes I\ = Wi from a function of p and the 
bare 0b to a function of p and the renormalized 0. This prescription fixes the counterterms 
and is analogous to that described for a 1PI effective action in Sect. 11.4 of Ref. [13J for 
tadpoles. Since the end result is simply the replacement 0b —* 0, there is no need to consider 
the explicit counterterms. We return below to the issue of additional finite contributions to 
( at this order. 

This renormalization prescription allows us to write W\ in Eq. f l4"6]) directly in terms of 
p and 0, yielding immediately the renormalized expression for I\: 



ri[p,0] = 7^V + 7^ 



(3V 



4 



With this result, we can evaluate pi and j\ directly: 



and 



Ji 



1 dVx 

/3V dp 



1 dT 1 



2°° 9 



(49) 



(50) 



(51) 



f3V 90 

More generally, we will need to eliminate p$ and jo from the Wi in favor of p and in order 
to construct the IV Since the particle number (and hence the density for a uniform system) 
is given, the correction pi is known. In contrast, the pairing density 0, which determines ji, 
is not fixed by Eq. fl38|) unless j (and p ) is known. But since we take j = in the ground 
state, and j = j + ji to this order, we have for the ground state that 



Jo 



-Ji 



2°° 



Substituting from Eq. (1381) . we get the "gap equation" 



Jo 



(i) 



3-D 



d D k 



o £ - C (0) Jo 



(27T) 13 £, 



a 



(1) 



-jo[/(o)-C 



(0)1 



(52) 



(53) 



where we've used 1(0) as defined in Appendix [B] This has the trivial solution j = but 
also the possibility of a non-zero "superconducting" solution for j Q (and therefore 0) for 
C«<0. 

If we evaluate the gap equation in DR/PDS using 1(0) from Eq. (|B5]) for D = 3 and C {0) 
from Eq. (140]) . the A dependence is precisely eliminated, leaving 



°o 



7(0) 



-^/2Mp7oWs\(l + x 2 ) 1/4 Pi/2 ( 



-1 



(54) 
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with IjO) defined in Eq. flB6]). With /x = kj/2M, this is the same gap equation derived in 
Ref. [40j, with j taking the place of the gap A. If we use the leading term in the expansion 
of 7(0) from Eq. (JB7J), we find 



x 



Jo 



Ik) 



2 C 



-it /2k-p \a s 



(55) 



as usual for the BCS weak-coupling limit. This is a very good approximation even for 
moderate values of fcp|a s |; in particular, the approximation for x is less than the exact 
solution to Eq. (PJ) by 0.06%, 0.8%, 2.7%, 3.8% for k F \a s \ = 0.5, 1.0, 2.0, 3.0, respectively. 
If we expand the expression for p in Eq. 
leading term cancels, leaving 



2J), we find that the logx dependence from the 



P 



h 3 

3^2 



1 + A x 2 (i_21og- 
16 V & 8 



Thus we get the free Fermi gas result with a correction that depends on x. For /cf|o s | k 
we find x ~ 0.05, and the correction is about 0.5%. 

The renormalized effective action to this order is found by combining Eqs. (I36p and 
and using the results in Appendix [B] 



(56) 
! 0.5, 



(3V 



(r„ + Ti 



2 



d D k 



(& - e„) + ^c (0) j 2 + pop - M + Ao'p 1 + -A 



(!) ^ 



1 



2 * -u ru, ,»t 4 

1AM 



- 7(2) + 2/i /(l) - + J 2 )n0) + j 277 
-7(2) + 2/i /(l) - 



■Jo + Pop - Jo</> 



(Po + Jo) AO) + Pop ~ Jo0 + \C { 1] P 2 + ±C< V 



(57) 



which is manifestly independent of A. To find the energy density, we evaluate at the sta- 
tionary point: 



'3o= 



J(2)+/i /(l) 



\a 



i 



-p 



\a 



(i) 



2 2 



(58) 



where we've used Eqs. (154]) and f|42l) to reach the final form. 

In practice, however, a different renormalization can be more useful for the numerical 
solution of the self-consistent equations and for extending the results to inhomogeneous finite 
systems in the local density approximation (LDA) and to higher orders. The conventional 
approach to renormalizing the gap equation for a delta-function interaction is based on the 
observation that the unrenormalized gap equation (with a gap A taking the place of jo in 
Ek), 



1 



I Co 



d 3 k 1 
WfE~ k 



(59) 



has the same linear ultraviolet divergence as the unrenormalized expression for zero-energy 
scattering in terms of the scattering length [39l . 51, 



M\C \ 

47TGU 



+ 1 



|C n 



d 3 k 1 
(2^F^ 



(60) 



17 



Thus, one can eliminate the unrenormalized Co directly (e.g., by subtracting the equations 
and dividing out Co) to obtain: 



M 


1 




f d 3 k 


" 1 


1 " 


4ira s 






(2tt) 3 







(61) 



which is finite and can be evaluated using elliptic integrals 1391. This same result is obtained 
within the present formalism by comparing Eqs. (|38|) and (I6T1) and making the association: 



(0) 



AM 



d 3 k 1 



(62) 



This connection gives us a way to eliminate explicitly in the energy and 4> equations, 
leaving finite expressions that can be evaluated numerically. 

Thus, the energy density from Eqs. fl5T|) and flo^l) can be evaluated as 



E 
V 



d 3 k 




-E k + 


lJo 2 l 


+ 


Po - 


4 l^o IP 


p + 


jI 


(2tt) 3 


2eJ 


l°o 1 


d 3 k 
(2vr) 3 


6 


-E k + 


IJo 2 
2E k 


+ 


Po - 


4 l^o IP 


p • 





(63) 



The expressions for p and in this approach are 

d 3 k 



P 



(27T) 



1 



■|r ) au<1 



-Jo 





" i r 


/ (27T) 3 


_E k Sk_ 



(64) 



We note that the usual normal-state LO results are recovered in the jo 

*M|„2 • 1 



limit: 



E 
V 



-Hop- ^\Co '\P 



and 



P 



37T : 



r A;| and fi 



kp 
2M 



(65) 



More efficient evaluations of have been discussed by Bulgac and Yu [55j, [56|, l57|] in the 
context of local density approximations (LDA's) for finite systems. In the present case, 
we simply illustrate how improvements can work by comparing the evaluation of for two 
different subtractions. We define <pE c as the evaluation of (f) as a numerical integral up to a 
maximum k c given by E c = k^/2M. Then we compare the two expressions for (pE c , 



-Jo 



r kc d 3 k 


" 1 1 1 




r kc d 3 k 


1 (2vr)3 


_Ek £k_ 


= ~ 30 J 


' (2tt) 3 



1 

Ek 



V 



£k — flQ, 



(66) 



where the subtractions are equal in the k c — > oo limit. However, as shown in Fig. [TJ the rate 
of convergence is, in practice, dramatically improved by using the second subtraction. This 
will be important in applications to finite systems. 

The simple interpretation of the effective action as a variational energy functional is 
apparently altered by the requirement of sources appearing beyond linear order In par- 
ticular, we recall the connection between an effective action and a variational estimate of the 
energy, either in Euclidean space [53[ or in Minkowski space as a constrained minimization 
with the sources acting as Lagrange multipliers [14| . We consider the former, with a generic 
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FIG. 7: Convergence of the integral for the pair density in the uniform system for two subtractions 
as a function of an energy cutoff E c in the integral. The converged value of <f> is denoted (p^. 



Hamiltonian H and external source J(x) with only spatial dependence coupled to a density 
f) (which represent vectors of sources and densities in the present case): 



H(J) = H + / Jp 



(67) 



If the ground state is isolated (and bounded from below), then as (3 — > oo the partition 
function in the presence of J, 2 [J], projects on the ground state of H(J) with energy 
Eo(J): 



Thus, separating out the explicit dependence on J in H(J), 



E (J) = Jim ~ log 2[J] = lim \w[J] . (68) 



E (J) = (H(J))j = (HI 



(69) 



where (O)j means the expectation value of O in the ground state in the presence of J. 
Combining these results, the expectation value of H in the ground state generated by J[p] 
is 



1 



(H)j = pW[J] 



Jf> = f?\p] 



En 



(70) 



Thus, we conclude that T[p] provides a variational estimate of the energy. However, it 
would seem that this schematic argument breaks down if the source J appears nonlinearly 
in Eq. (ISTl) . Similarly, the interpretation of J as a generalized Lagrange multiplier relies on 
it appearing linearly [14J ] . Another alternative would be to remove nonlinear terms in favor 
of auxiliary fields, as discussed in Refs. (H, [13]. However, it is still the case that the ground 
state energy is given at a stationary point (corresponding to J — > 0), as long as the operator 
corresponding to p is defined through p = 6W[J]/5J. 
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V. RENORMALIZATION TO NEXT-TO-LEADING ORDER 



In this section, we extend the inversion procedure to next-to-leading order (NLO) at zero 
temperature, which means calculating r 2 [/io, jo] an d determining equations for fi and j . 
We start with the diagrams for W% in Fig. []Jb) and (c) evaluated with Cq vertices, together 

(2) 

with the W\ diagram evaluated with the C Q vertex. The latter diagram has two terms, each 
of which will remove a linear A dependence corresponding to ultraviolet divergences from 
each of the two W2 diagrams. The first term cancels a divergence and A dependence from the 
beachball, just as in Sect. Ill CI The second term cancels a new divergence in the anomalous 
diagram. Finally, we will have the counterterm proportional to j(iptip\ + which is 

again chosen to convert 0b to the renormalized in the anomalous diagram Fig. QJc). 

Applying the Feynman rules to W 2 an d and then carrying out the frequency 

integrals and simplifying through variable transformations, we obtain 

1 :W?WM - --Ac^ 2 [ dip fd4k fdS 



py-z Lr-u,,^ 4 w ) J (2?r) 4 y (2 7T) 4 7(27r) 4 

x {Tr[G (p)r 3 Go(p - q)r 3 ] ■ Tr[G (k)r 3 G (% + q)r 3 
- Tr[G (p)r 3 G (p - q)T 3 Go(k)r 3 G (k + q)r 3 ] 

+ Top T lp-q Ftf^T 3k + q + 2 Gop Qop-q ^0% ? 3fc+g 

d 3 p f d 3 k f d s q 1 



•(ci 



(2tt) 3 J (2tt) 3 J (2tt) 3 E p + E k + E p _ q + E k+q ' 
x [u 2 p ulvl_ q vl +q - 2n\v\ (uv)p-q (uv) k+q 

+ (uv) p (uv) k {uv) p - q (uv) k+q ] (71) 



and 



:WfWo] = l(C^ fdP fdk fdq 



py* WW 2 v-o j jj^yL JJ^y iO^r) 1 

x {Tr[G (p)r 3 ] ■ Tr[G (fc)r 3 G (fc)r 3 ] • Tr[G (g)r 3 ] 

- 2Tr[G (p)r 3 ] • Tr[G (k)r 3 G (k)r 3 G (q}T 3 } 
+ Tr[Go(p)r 3 Go(fc)r 3 Go(g)r 3 G (fc)r 3 ] 

d 3 k 1 r , x9 1 



-(C«)^^^[pK^) 2 + ^bK-^)] 2 ■ (72) 



We can identify potential ultraviolet divergences by noting that 



«-£(!- (73) 



4 



2( 1 + ^)- J?1 -V (74) 



20 




FIG. 8: Renormalization of ultraviolet divergences at NLO. The divergent parts of W2 evaluated 
with Cq 1 ^ are shown with the corresponding W\ contributions evaluated with , which precisely 
cancels the A dependence. 



jo fc->oo j M 

UkVk = -2E k ^-^ (75) 

and 

1 k ^oo2M 

E k - IF ■ (76) 

Using these asymptotic behaviors, we see that in W 2 it is only the term with u\ui that 
can be ultraviolet divergent (see Fig. [Sfe.), while in W% it is the term with u\ times 0| (see 

fig. Eb). 

To isolate and cancel the divergent contribution in W^, it is easiest to switch variables 
to s, t, and u, where p = s + u, k = s — u, p — q = s — t, k + q = s + t, and we pick up a 
factor of 8 from the Jacobian. After replacing ul± u by 1 — v 2 ± u and keeping only the factors 
with l's, we find the divergent part 

(3V [W2 Jdlvergent " 8(C ° } J (2tt)» J (2vr)3 J (2vr)3 E s+U + E s _ u + E s+t + E s _ t " ^ 

For large w, we see from Eq. (1761) that Eq. ( 1771) has the same linear divergence as found for 
82 in the jo — limit in Sect. Ill CI We stress that with the divergence regulated in DR, 
consistent renormalization requires that the A dependence be removed by the contribution 
from the Hartree-Fock diagram evaluated with the C vertex (but now with j 7^ 0) without 
adjustment. The latter is found from the p 2 term in W\ [see Eq. (H6|) ] with — > = 
47ra 2 s A/M, 

1 4Tra 2 s A f dH [ d 3 k> 2 , 2 

w m[ * = ^JwJw** ' ■ (78) 

To carry out the same sort of DR analysis as for LO, we would need to analytically continue 
Eq. (177}) in D. A more practical alternative is to subtract and add the jo = version of the 
integral (also taking e s and St less than po), which can be evaluated using the integral in 
Eq. (US]): 

-^ 1} ) 2 /A /A (ST"* f d ° u 



(2tt) 3 J (2tt) 3 V 2 ; j (27r) fl e s+u + e s - u - e s+t - £ s _ t 

,(l)x2 / d 3 S f dH 2 2 {A\ 3 ~ D /" gPw P 



47ra?A /" <i 3 s f dH 9 9 , . 

^s-t^s+t, (79) 



M j (2vr) 3 y (2vr) 3 
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which, after a change of variables, is seen to equal — ^<W 2 . Thus, the renormalized version 
of Eq. (1771) is obtained with the substitution: 

1 1 V 



E s + U + -^s-u + Es+t + ^s-t E s+U + E s _ u + -E s+ t + E s _ t £ s +u + ^s-u ~~ ^s+t ~~ £s-t 

(80) 

which makes the integral explicitly finite and suitable for numerical evaluation. 

(2) 

The other part of W% with C§ renormalizes the anomalous W 2 contribution. In par- 
ticular, we isolate the divergence in Eq. (172]) by the replacement u\ — > 1 — v k . Then the 
divergent part is 

>*V- = -(^)7(g 5 ^^- (81) 
In this case, the substitution 

1 1 ^ ( 82 ) 

•E* E k e k 

renormalizes the anomalous W 2 . Using Eq. (jS2"|) . we see that the additional piece is supplied 
precisely by the second term in W\ with Cq — > Cq 2 \ as advertised. 

The prescription for deriving T 2 from W 2 was detailed in Ref. [9j for the case of a source 
J(x) coupled to the density, and is reproduced for a constant chemical potential in Eq. ( fLfl) . 
The end result is simple: the additional terms from the inversion method cancel precisely 
against the anomalous diagram (for zero-range forces). This same cancellation holds more 
generally with multiple sources and corresponding Legendre transformations (the demon- 
stration follows directly by putting sources and densities into a vector, and is not shown 
here). Here one has to check that the counterterm contributions are consistent with this 
cancellation. After the removal of the anomalous diagram W 2 , the evaluation of T 2 [p, 4>] 
reduces to evaluating the renormalized beachball diagram. That is, T 2 is given by Eq. flTTj) 
with the term corresponding to Eq. (1771) evaluated with the substitution in Eq. f[8"Uj) . 

While the NLO result for the energy follows formally as at LO, in practice the solution 
does not follow as directly. It particular, analytic expressions are not available for the T 2 
integrals and j 2 = — (1/ j3V)dT 2 /d(p\ p is not simply proportional to (unlike ji), and so 
jo = —ji — j 2 does not yield a simple gap equation. A possible solution method is to directly 
solve <9[r + Ti + r 2 ]/d(j)\p = 0, for example via a steepest descent method, by varying /i 
and jo while keeping p constant using Eq. fPf2l . This task can be simplified by expanding 
about the LO solution and/or using the expansions for small x = jo/po- Such numerical 
solutions have yet to be investigated. 

Although the full result for T 2 [p, 0] is only accessible numerically, we can extract the 
leading correction to the gap equation analytically, which at NLO results in a universal (i.e., 
independent of a s ) change in the prefactor of jo (which we associate with the gap). This 
contribution can be extracted by isolating the contribution to T 2 that has the logarithmic 
divergence at the Fermi surface of the same form as in W\, which means two factors of uv. 
We can identify this part in W 2 ^ by letting p — q — > k and k + q — > k' in T 2 to obtain 

— [^Jinduced = 2(6^) J7T-^U k V k I J7rT^ u k' v k' 1 



(3V L " u v u 'J (2tt)3 * V (2vr)3 * K J (2tt)3 £ k + £ k , + £ k+q + £ k ,^ q 

(83) 

We note that the terms proportional to u k v k peak at p§ as jo — > (see Fig. [9]). To extract 
the leading contribution from the log divergence at the Fermi surface, we can set jo = 
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FIG. 9: Left: Plot of u\ and v\ along with UkVk- Note: the width labeled "jo" is exaggerated to 
make it more visible. Right: Same plot with k 2 v\ (dotdashes) and k 2 UkVk added, with an extended 
x-axis. 



with |k| = |k'| = k-p in the inner integral. This means that 27 k , i? k / — > and the u 2 and v 2 
terms become theta functions. As suggested by Fig. [TOJ, the inner integral has the form of a 
non-interacting density-density correlator or Lindhard function. Such a correlator is defined 
and evaluated in Ref. [45[ for non- interacting T = Green's functions Gq\ 



n°(q,g ) = -2i 
= 2 



d A k 



G°(k)G°(k + q) 



{2*y 

d 3 k 

^6(\ q +k\-k F )6(k F -\k\) 



x 



1 



,go + £k-£ q +k + ^ ?o + £q+k-£k-^ 
Because of the peaking at /io, the leading contribution is 



[-^2] induced 



l/^,fiis9 f d 3 k f d 3 k' „ n/ , , , , 1 



(84) 



|k|=|k'|=fc F 



(cf )VI (n (k-k'))| kNk , 



|k| = |k'|=fc F > 



(85) 



where (■ • •) denotes an average over the angle between k and k'. Using the expressions in 
Ref. [45j and defining 77 = a/ (k F — k ■ k.')/2k F , 



(n (k - k'))| k | = | k /| =fcp 



Mfcp 1 
~2^~3' 



drj 



' , 1-77% /! 
rj H in 



1 + 2 In 2) 



From Eq. (1551) . the leading order (LO) result for |jo/)Uo| at T = is 



(86) 



Ij'o/m 



LO — ~R 

e z 



-7r/2fcp|a s | 



(87) 
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ri + r 2 = * +41 + 




k u k 

FIG. 10: Contributions to Ti and T2 with logarithmic divergences at the Fermi surface as jo — > 0. 
The vertices in each case have coefficient Cq 1 ^. 

After adding the renormalizations to T 2 that take </>b to <f), we find that Eq. fl52l) becomes 
at NLO 

jo = -(h+h) = -h c o ] \ 



\ , c (i),MAft / l + 21n2 



2tt 2 



(88) 



The final result for |j'o/a*o| is that the exponent is modified, which changes the prefactor, 
yielding 

Uo/HoInlo ~ (4 e )i/3 Uo/^oUo ( 89 ) 
or, associating jo with the gap again (see below about the gap in a Kohn-Sham DFT), 

A NLO « ^TTs A ^o • (90) 



Thus, we recover the universal suppression of the gap in the kp — > limit found in Refs. [41 
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VI. SUMMARY 

In this paper, we extend the effective field theory treatment of a uniform, dilute Fermi 
gas to include pairing. The starting point is a generating functional for the grand canonical 
partition function including a chemical potential for the fermion number and a source coupled 
to the pair density, which we take to be constant. Legendre transformations with respect to 
the sources lead to an effective action of the fermion and pair densities. 

The effective action formalism is often used in condensed matter applications to discuss 
superconductivity. In that case, the pairing field enters as an auxiliary field. After integrat- 
ing out the fermion fields, a conventional one-particle-irreducible effective action is derived 
by Legendre transformation with respect to the auxiliary field. A minimum at nonzero 
expectation value of this field (which is proportional to the pairing gap) indicates sponta- 
neously symmetry breaking of the phase symmetry related to fermion number conservation, 
which implies the normal ground state is unstable to pairing. 

In contrast, the calculation here is carried out by adapting the inversion method pro- 
posed in Ref. pi| to an EFT treatment of the dilute Fermi gas. The inversion method in 
Ref. [l9| was applied to conventional BCS superconductivity using a source coupled to the 
pair creation and destruction operator. Instead of organizing the perturbative inversion in 
terms of a coupling constant (e.g., the electron charge squared), we use the EFT expansion 
parameter for the natural dilute system, which is the inverse of the resolution scale. In the 
uniform system, the ultimate dimensionless expansion parameters are products of the Fermi 
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momentum kp and parameters of the effective range expansion (e.g., the scattering length 
a s , which is the order of the interaction range for a natural system). 

We encounter a new renormalization problem because of ultraviolet divergences involving 
the composite pair density. To deal with it, a new term proportional to the square of 
the external source j, with new coupling constant (, is added to the Lagrangian. In our 
renormalization prescription, we implicitly assume that ( has a Taylor expansion in a s to be 
consistent with the inversion method. On dimensional grounds, a term with ( proportional 
to l/a s is possible, which would be independent of A. Such a term is consistent with the 
development for effective actions with local composite operators in Refs. jH, SI, where 



the coefficient of analogous terms have a Laurent expansion starting at the inverse of the 
coupling. However, the counterterms in our prescription along with the source-pair-density 
counterterms suffice for us to renormalize, with the prescription that the pair density is 
unchanged from the non-interacting system. We verify the renormalization through NLO 
by using the PDS subtraction prescription with dimensional regularization. 

At leading order, we reproduce the results of Papenbrock and Bertsch [!o[ for the gap 
and energy density. At next-to-leading order, we reproduce the leading weak-coupling limit 



derived by Gor'kov and Melik-Barkhudarov [41( for the modification of the gap, and generate 



complete expressions for the observables suitable for numerical evaluation and application 
to finite systems in the local density approximation (LDA). We have focused on the zero 
temperature limit here, but the formalism applies directly at finite temperature. For exam- 
ple, we can find the transition temperature by looking for where the non-trivial (jo 7^ 0) 
solution to the j self-consistency equation disappears, which corresponds to solving for the 
temperature with j set to zero. 

The next step is to extend our results to finite systems. This means merging the formalism 
here with the density functional formalism of Refs. 0, [Tol. [Tlj] . The renormalization will carry 
over from the present discussion via the LDA, as in Ref. [9j. In practice, efficient subtraction 



schemes at leading order have been developed by Bulgac and Yu |55|, |56|, |57[ . 

The extension to higher-order terms, such as C2 and C 2 requires further renormalization. 
In general, we expect that every term of equal or lower dimension, consistent with sym- 
metries, to enter. When C2 is added, for example, in the general case a new source term 
proportional to C^j 2 "?/^ is needed. Extending to isospin and spin-dependent interactions 
is straightforward formally, but the renormalization issues should be reconsidered at each 
stage. 

We have noted that the Kohn-Sham potential jo plays the same role at leading order that 
a gap does in a BCS treatment with a contact interaction. The conventional definition of 
a gap function is an integral over the potential and the anomalous Green's function. For 
example, Mahan defines for a local potential V^(q) [54j : 

A(q) = "7 E ^( k )^o(q - k, r = 0) = - J 0^ V(q - k)F (k , k) e *" . (91) 

Given a momentum expansion of V, we could define a corresponding series of gap functions 



(see the appendix of Ref. 40]). However, it is important to note that the Kohn-Sham "gaps" 
are not directly connected to the single-particle spectrum. In Ref. [Tlj], the connection 
between the Kohn-Sham Green's function and the full single-particle Green's function is 
discussed. There is no difference at LO, but at NLO and higher orders the self-energies will 
differ, which means there is no longer a direct correspondence between the Kohn-Sham gap 
function and that extracted from the full Green's function. The systematics of the differences 
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between these gaps (and also a gap denned from even-odd staggering of ground-state energies 
in finite systems) have yet to be investigated. 

The present formalism is not the only path to a nuclear DFT including pairing. As 
discussed in the introduction, an alternative to a local source coupled to the local pair 
density is to use a nonlocal source and density 0, 29[. Proceeding in this way would 



parallel the response in particle theory to problems uncovered by Banks and Raby, which 
was to focus on CJT two-particle-irreducible effective actions. An alternative would be to 
apply the auxiliary field approach in the pairing channel. Here a clean separation of particle- 
hole and particle-particle channels may be the problem, although this is regularly a part of 
phenomenological approaches. Finally, a completely different path to DFT is outlined by 
Polonyi and Schwenk, who propose a renormalization group method [58]. These avenues are 
all worth exploring. 
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APPENDIX A: CANONICAL DERIVATION OF W [// , jo] 

In this appendix, we derive an expression for W [fi , j ]/ j3, which is proportional to the 
logarithm of the partition function of the noninteracting system in the presence of /xo and 
jo, in the zero temperature {(3 — > 0) limit. If /x ^ but j = 0, we project onto the 
normal noninteracting ground state, which is simple to construct. In the uniform case, we 
just fill the lowest plane- wave states (which have the lowest eigenvalues) up to the Fermi 
momentum kp, which is related to the density by Eq. (TTTT) . This minimizes the energy for 
a given particle number. When jo 7^ 0, we need to construct the state that minimizes the 
free energy for given chemical potential /xq (and given jo). Here we follow the notation and 



general discussion of Ref. [45] in using the Bogoliubov transformation, which finds the exact 
diagonalization in the noninteracting case with constant jo- Not surprisingly, the state that 
minimizes K Q = H — /xo-/V, where H includes the jo dependence, is a BCS-type variational 
state. 

The non-interacting Kq in terms of creation and destruction operators for the momentum 
states with spin up and spin down is: 

^0 = 5>k - Po) «L«k Q + X> (4 T «-k| + «-k|«k T ) + ^C (0) J 2 V , (Al) 

ko k 

with = k 2 /2M. (This definition generalizes K = H — fi N from Ref. (45[.) This result 
is apparent from inspection, but can be derived formally from the field expansions for ^(x) 
and V^( x )- The term proportional to jg is needed to ensure A-independent results, but plays 
no explicit role here. The minimized expectation value of K at fixed V and /x at T = is 
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related to the energy by 

^WoKio] ^ (K ) = n ( f i,T = 0, V) = (E - HoN)\ T=0 . (A2) 

We introduce the Bogoliubov-transformed operators a k and /3 k , 

«k = MfcO kT - v k a ] _ n , /3 k = M fe a„ ki + f fe a kT , (A3) 

so that 

a kT = MfcCtk + v k pl k , a_ kj = w fc/ 3_ k - v k a[ . (A4) 

The anticommutation relations 

{a k ,a k ,} = {/3 k ,/3 k ,} = (5 kk , , (A5) 
with all others equal to zero, imply that 

ul + v 2 k = l. (A6) 
The BCS ground state is defined so that 

a k |BCS) = A<|BCS) = . (A7) 



(See Ref. 45( for an explicit construction in terms of the a's and cr's.) At this point we 
don't know that this is actually the ground state of the noninteracting system with non-zero 
jo, but we do not need to make this assumption. 

By following the normal-ordering approach in Ref. jiBJ or simply by applying the com- 
mutation relations, we can rewrite K with no loss of generality as: 

K = U + H 1 + H 2 , (A8) 

where 

U 



2 Zk4 + 2jo UkVk + V » ( A9 ) 

k k 

Hi = ^(«k«k + pt v P-*)[(ul - v 2 k )Ck - 2u k v k j ) , (A10) 

k 

H 2 = X)(4^k + P-*<Xk)[2u k v k Z k +j Q (ul - v 2 k )) . (All) 

k 

In these expressions, 

& = el-Lio. (A12) 

We still have freedom in our choice of u k and v k , subject to the constraint (IA6I) . Using this 
freedom, we can choose 

2u k v k £ k = -{u\ - vl)j , (A13) 
so that H 2 vanishes identically. This condition is satisfied for 

^=-4- --K i+ I)' (A14) 
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with 

Ek = \/e k +fo- (A15) 



Substituting from (1A14I) into Eq. flAlOp . we find that Hi simplifies to 

h x = ^*(4«k + /3i k /u) • (M6) 

k 

Then it is clear from this equation together with Eq. (IA70 that the corresponding |BCS) 
state (which is dependent on /i and jo) will be the ground state of A" (/i , j ), since U is a 
c-number, Ef. is positive definite, and (H 2 ) = 0. 

Note that there is no self-consistent gap equation at this point; for any given value of /xo 
and jo, w e can construct this ground state with Eq. (IA14I) satisfied. The density of fermions 
follows immediately as 

"=y = ^E< BCS l4AjBCS) = |i:^=Ii:(l-|) . (A17) 

ka k k V kJ 

and the ground-state expectation value of K , which is equal to U, is 



U = (BCS|A^o|BCS) = I ^W [fx , Jo 



P 

2 tkvl + 2j'o u ^ + ^ (0) io V 



= E^-^) + ^ (0) jo 2 ^. (Ais) 

k 

Upon converting the sums to integrals and dividing by the volume, we obtain the expression 
for Wq/ (3V in Eq. (1331) . We emphasize that there is no approximation or truncation here. 

APPENDIX B: INTEGRALS FOR DR/PDS 



Here we generalize the DR/MS integrals used in Ref. |40j to the PDS subtraction scheme. 
The starting point is the integral [4(| 



7( = 2 = (a 2 + br /2 (—) P° a (-a/Va^) , (Bl) 

Jo ^{e-af + b 2 \sm7raj 

where P° is a Legendre function. The right side of Eq. ( 1B1I) can be analytically continued in 
a to define the relevant integrals in D dimensions. For appropriate D, there is an analytic 
formula, and then the answer is continued back to the relevant number of spatial dimensions 
after the approprite subtraction of poles (according to the prescription). 
A particular integral of interest with a = is 

/ A\ 3 ~ D f d D k 1 /2AA D/2 (A/2) 3 - D r°° e^/ 2 " 1 ) de 

W J (2vr)V^k-/io) 2 Tl ~ P2) Jo T^oT+l 
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^2Mfi \ D/2 (A/2) 



47T 



3-D 



-7T 



r(£>/2) Vsin^l^ 1 / 2 - 1) 



x 



(l + x 2 )( D/2-l)/2 



Mo 



D/2-1 



(B2) 



where x = Ij'o/MoI- To apply Eq. (1B2[) in the PDS scheme, we must subtract poles in (2 — D). 
The only such pole comes from 



-7T 



7T 



+ — (2 - Z>) + 0[(2 - D) 



sin7r( J D/2 - 1) 2 - D 12 
which means the pole term is (all other factors evaluated at D = 2): 

2 \ 2M/i MA 1 D-> 3 MA 



-A 
IT 



1 

Mo 



2-DJ 4vr 2tt (£>- 2) 2vr 

More generally, the basic DR/PDS integral in Z) dimensions, with x = Jo/Mo; is 

MA 



A X 3-z? /• d D k (el) 



+ 



(2vr)^ £ fc 
0+1 ^ 3/2 



(Mo)' 3 ( 1 



V2- 



7T 



2vr w V " 2 

[^(l + x^f+WpO 



/3+1/2 



1 



MA 

~2tT 



(Moffl-^y )+/(/?) 



where 7(/3) is the A-independent part of I(/3) 

We can expand /(/?) with f3 = 0, 1, 2 for small x as: 



7(0) 

7(2) 



(2M/i ) 1/2 M 



7T^ 



x 



2 + log-) +^ (l + log| N ) 



(2M/i ( 



\3/2 



2tt 2 
(2^o) 3/ Vo 



+ log 



.r 



log 



X 



2tt 2 



46 , x 
h log — 

15 S 8 
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x 2 + 



X 



1 + log - )x z + 



(B3) 



(B4) 



(B5) 
(B6) 
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